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ABSTRACT 

Context. Imaging of debris disks has found evidence for both eccentric and offset disks. One hypothesis is that they provide evidence 
for massive perturbers, for example, planets or binary companions, which sculpt the observed structures. One such disk was recently 
observed in the far-IR by the Herschel Space Observatory around // 2 Reticuli. In contrast with previously reported systems, the disk 
is significantly eccentric, and the system is several Gyr old. 

Aims. We aim to investigate the long-term evolution of eccentric structures in debris disks caused by a perturber on an eccentric orbit 
around the star. We hypothesise that the observed eccentric disk around f 2 Reticuli might be evidence of such a scenario. If so, we are 
able to constrain the mass and orbit of a potential perturber, either a giant planet or a binary companion. 

Methods. Analytical techniques were used to predict the effects of a perturber on a debris disk. Numerical (V-body simulations were 
used to verify these results and further investigate the observable structures that may be produced by eccentric perturbers. The long- 
term evolution of the disk geometry was examined, with particular application to the / 2 Reticuli system. In addition, synthetic images 
of the disk were produced for direct comparison with Herschel observations. 

Results. We show that an eccentric companion can produce both the observed offsets and eccentric disks. These effects are not imme- 
diate, and we characterise the timescale required for the disk to develop to an eccentric state (and any spirals to vanish). For f 2 Reticuli, 
we derive limits on the mass and orbit of the companion required to produce the observations. Synthetic images show that the pattern 
observed around ( 1 Reticuli can be produced by an eccentric disk seen close to edge-on, and allow us to bring additional constraints 
on the disk parameters of our model (disk flux and extent). 

Conclusions. We conclude that eccentric planets or stellar companions can induce long-lived eccentric structures in debris disks. 
Observations of such eccentric structures thus provide potential evidence of the presence of such a companion in a planetary system. 
We considered the specific example of ( 2 Reticuli, whose observed eccentric disk can be explained by a distant companion (at tens 
of AU) on an eccentric orbit (e p > 0.3). 

Key words, circumstellar matter - methods: numerical - planetary systems - stars: individual: /" 2 Reticuli - celestial mechanics 


1. Introduction 

The first debris disk was discovered in 1984, when the 
Infrared Astronomical Satellite (IRAS) found a strong IR excess 

* Appendices are available in electronic form at 
http : //www . aanda . org 

** Herschel Space Observatory is an ESA space observatory with sci- 
ence instruments provided by European-led Principal Investigator con- 
sortia and with important participation from NASA. 


around Vega, revealing the presence of micron-sized dust grains 
(Aumann et al. 1984). For most debris disks, these grains have a 
limited lifetime, which is shorter than the system’s age because 
of Poynting Robertson drag and collisions. Therefore, this dust is 
assumed to be replenished by collisional grinding of much larger 
parent bodies, which are at least kilometre-sized for this colli- 
sional cascade to be sustained over the system’s age (Backman 
& Paresce 1993; Lohne et al. 2008). Consequently, debris disks 
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provide evidence for the existence of solid bodies that have 
reached km-size, and potentially the planetary-size level. 

Spatially resolved structures in debris disks can provide 
clues to the invisible planetary component of those systems. 
Such planets may be responsible for sculpting these disks and 
may leave their signature through various asymmetries such 
as wing asymmetries, resonant clumpy structures, warps, spi- 
rals, gaps, or eccentric ring structures (see, e.g., Wyatt 1999). 
This diversity is to be compared with the variety of exoplane- 
tary systems 1 discovered around main sequence stars since 1995 
(51 Peg b. Mayor & Queloz 1995). Dynamical modelling of 
such asymmetries is the only method to place constraints on the 
masses and orbital parameters of planets in systems where di- 
rect observations are not possible (see, e.g., Mouillet et al. 1997; 
Wyatt et al. 1999, 2004; Augereau et al. 2001; Moro-Martin 
& Malhotra 2002; Kalas et al. 2005; Quillen 2006; Stark & 
Kuchner 2008; Chiang et al. 2009; Ertel et al. 2011; Boley et al. 
2012; Ertel et al. 2012; Thebault et al. 2012). 

We focus here on cases of eccentric patterns in debris disks. 
The modelling of this type of asymmetry and its possible link 
with the dynamical influence of eccentric companions has been 
investigated in several previous studies: authoritative work was 
carried out by Wyatt et al. (1999, 2000) for HR 4796. Another 
case of interest is the debris disk of Fomalhaut (Stapelfeldt et al. 
2004; Kalas et al. 2005, 2013; Quillen 2006; Chiang et al. 2009; 
Boley et al. 2012; Kalas et al. 2013; Beust et al. 2014). 

This pioneering work showed that these large-scale struc- 
tures arise in systems where debris disks are perturbed by an ec- 
centric companion on a low-inclination orbit relative to the disk 
(Wyatt et al. 1999). The disk centre of symmetry is offset from 
the star, which may be measured explicitly in high-resolution 
images (e.g., HST scattered light). Furthermore, its periastron is 
closer to the star and thus hotter and brighter, which results in a 
two-sided brightness asymmetry. 

However, it is important to emphasize that previous 
studies considered low-eccentricity rings ( e > 0.02 for 

HR4796 and e = 0.11 ± 0.01 for Fomalhaut), and were lim- 
ited to timescales shorter than the typical ages of mature 
disks (<10 Myr for HR4796 simulations and <100 Myr for 
Fomalhaut). The question of whether highly eccentric ring struc- 
tures could be sustained over very long timescales has not 
been addressed thus far in the literature. This question has be- 
come very topical because of the discovery of at least two 
several Gyr old and significantly eccentric debris disks: one 
around f 2 Reticuli ( e > 0.3, Eiroa et al. 2010), which is used as a 
reference case in this paper, and another one around HD 202628 
(e ~ 0.18, Stapelfeldt et al. 2012; Krist et al. 2012). These sys- 
tems are older than either Fomalhaut or HR 4796, and their disks 
are also much more eccentric. 

In the present work, we investigate the long-term evolution 
of highly eccentric structures in debris disks and their relation 
to planetary or stellar perturbers by investigating their evolution 
over Gyr timescales. One of the questions we address is whether 
these structures are indeed several Gyr old, or if they might have 
originated from a more recent event, either a flyby or the late 
excitement of a sherpherding planet’s eccentricity. We also sum- 
marise a general modelling method, based on complementary 
analytical and numerical tools, which we apply to the specific 
case of l 2 Reticuli. 

This paper is structured as follows: Sect. 2 presents how a 
perturber can generate an eccentric ring structure. Useful ana- 
lytical expressions are derived to study under which conditions 

1 See http: //www. exoplanets . org 


such a pattern can be created. We also show that these predic- 
tions can be complemented by numerical studies. Section 3 de- 
scribes the debris disk of £ 2 Reticuli, along with newly reduced 
HerschelfPACS images. This debris disk is used as a proxy to 
determine a numerical set-up. Then, in Sect. 4, the numerical 
investigation is carried out. From A-body simulations, we exam- 
ine both the onset and suvival of an eccentric pattern and explore 
their dependencies on the perturber’s characteristics. This mod- 
elling approach allows one to put constraints on a perturber at 
work in shaping a debris disk into an eccentric ring over Gyr 
timescales, and it is applied to the debris disk of C 2 Reticuli. 
Section 5 shows synthetic images on which we perform a full 
comparison with observations of f 2 Reticuli, and retrieve ad- 
ditional constraints on this disk. Finally, Sect. 6 is devoted to 
conclusions, discussions, and propositions for future work. 

2. Footprints of eccentric companions on debris 
disks 

We have developed a dynamical model to investigate the shaping 
of a debris disk into an eccentric ring, and the timescales associ- 
ated with its onset and survival. More specifically, we seek to de- 
termine whether perturbers are able to shape and maintain a disk 
into a significantly eccentric ring structure on Gyr timescales, 
and whether the asymmetry relaxes or not. 

Before presenting our model and our results in detail for this 
as yet unexplored case, we present the current understanding on 
how eccentric ring structures arise as a result of the dynamical 
effect of an eccentric perturber. 

2 . 1 . Basic principle 

For a disk to be shaped into an eccentric ring, it must be per- 
turbed in such a way that the eccentricity of its components are 
forced to higher values, and that the orbits of these components 
are more or less oriented in a common direction. These con- 
ditions are both fulfilled if the disk is under the gravitational 
influence of a perturber, namely a planetary or a stellar com- 
panion (nearly) coplanar to the disk and on an eccentric orbit. 
The eccentricity of the ring causes the disk centre of symmetry 
to be offset from the star, and the disk pericentre to be brighter 
than the apocentre, since it is closer to the star and thus hotter. 
This feature was studied by Wyatt et al. (1999) for a planetary 
companion, and was called the pericentre glow phenomenon. 

Spatially resolved imaging is required to determine the struc- 
ture of debris disks, and therefore renewed efforts have been 
made to image as many debris disks as possible 2 . Most images 
of resolved debris disks have been obtained so far in the visi- 
ble or near-IR. At these wavelengths, the emission is dominated 
by small grains close to the blow-out limit imposed by stellar 
radiation (sub-micron to micron, depending on stellar luminos- 
ity). These grains exhibit a complex evolution because of the 
coupled effects of collisions and radiation pressure (see, e.g., 
Thebault & Augereau 2007). This may strongly alter or even 
mask the dynamical structures imparted by a massive perturber. 
Observations at longer wavelengths detect larger grains that are 
less affected by stellar radiation (a few tens to a few hundreds of 
micron in size, depending on the observing wavelength). These 
directly trace the distribution of larger parent bodies and thus 
more directly reflect the dynamical effect of a perturber (see, 
e.g., Krivov 2010; Moro-Martfn 2013, for exhaustive reviews). 


2 See, e.g., http://www.circumstellardisks.org 
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Since the origin of an eccentric pattern is gravitational, we 
can reasonably assume that large-scale asymmetries among an 
observed dust population already exist amongst the parent plan- 
etesimal population that produces it and result from pure grav- 
itational perturbations. This assumption allows us to study the 
influence of different eccentric perturbers in a simplified way, 
neglecting the effect of radiation pressure and considering par- 
ent planetesimals as mass-less and collision-less particles in or- 
bit around their host star and perturbed by a companion, either 
stellar or planetary. 

We assumed that at the end of the protoplanetary phase, the 
planetesimals start from almost circular orbits because of or- 
bital eccentricity-damping by primordial gas. We furthermore 
assumed that any perturbing planet in the system is fully formed 
by the time the gas disappears, and evolves on a significantly ec- 
centric orbit, because of a major perturbing event such as planet- 
planet scattering. Thus, we can consider the disappearance of 
the gas as time zero for the onset of planetesimal perturbations 
by an eccentric companion. From this moment, the planetesimal 
eccentricities start to increase and their lines of apsides tend to 
align with this of the planet. Under these assumptions, the forced 
elliptic ring structure takes some time to settle in, and it is pre- 
ceded by the appearance and disappearance of transient spiral 
features. These are due to differential precession within the disk: 
all the planetesimals in the disk have different precession rates 
(because of their different orbital distances), such that these spi- 
ral structures are expected to wind up and finally generate an 
eccentric ring, as shown by Augereau & Papaloizou (2004) and 
Wyatt (2005). The characteristic time for reaching this state is of 
the order of a few precession timescales at the considered dis- 
tance (Wyatt 2005). Consequently, the onset of an eccentric ring 
structure is a matter of timescales, while the value of the disk 
global eccentricity is to be linked with the planetesimals’ forced 
eccentricity, and thus to the companion’s eccentricity. 

2.2. Analytical approach 

We show here how the onset of these structures can be under- 
stood from analytical considerations. Planetesimals are consid- 
ered to be mass-less particles. We focus on the secular response 
of a debris disk to a coplanar perturbing body, either a planet or 
a star. More specifically, both the forced secular eccentricity et 
and the apsidal precession rate dm/dt of test planetesimals are 
examined, where et is the longitude of periastron with respect 
to the direction of the perturber’s periastron, that is, the planet 
and planetesimal have their periastra aligned when et = 0 and 
anti-aligned when et = n. 

When secularly perturbed, the eccentricity of a planetesi- 
mal evolves cyclically; its period is related to the rate of or- 
bital precession. This holds in particular when we consider a 
dynamically cold disk of planetesimals as an initial condition, 
which is a reasonable and classical assumption considering the 
damping effect of the gas during the protoplanetary phase. In 
that case, the secular behaviour of a planetesimal can be under- 
stood considering the analytical solution for its eccentricity. In 
the Laplace-Lagrange theory, the complex eccentricity of a plan- 
etesimal, z(t), can be written 

z(t) = e f {1 -exp(Mf)}, (1) 

where I 2 = - 1, and A = ^ is the secular precession rate. 

One can see from this expression that the maximum induced 
eccentricity for a planetesimal is twice the forced eccentricity, 
that is, ef >max = 2ef. This occurs when At = n[2n\, that is, when 
the longitudes of periastra of the planetesimal and the perturber 



AU 


Fig. 1. Co-evolution of a planetesimal eccentricity and orbital preces- 
sion when acted upon by an eccentric perturber. Since the planetesi- 
mal orbit precesses, its periastron will eventually be aligned with that 
of the perturber on the same side of the massive central body, that is, 
the longitudes of periastra of the planetesimal and the perturber are 
equal (et = 0). In this configuration, the planetesimal eccentricity is 
maximum. 


are equal (et = 0, see Fig. 1 and, for more details, see e.g., Wyatt 
2005; Beust et al. 2014). 

There are several ways to analytically derive ef and dm/dt. 
The most classical one is to apply the linear Laplace-Lagrange 
theory, which is an expansion of the interaction Hamiltonian to 
second order in ascending powers of the eccentricities of the two 
bodies and an averaging over the two orbits (see, e.g., Eq. (6) of 
Mustill & Wyatt 2009). However, this approach is valid only for 
low eccentricities, but the perturber’s orbital eccentricity e p is 
not necessarily low. Therefore, restricting our analytical study to 
small e p may not be appropriate. 

Another way to proceed is to expand the interaction 
Hamiltonian in spherical harmonics and truncate it at some or- 
der in a, where a is the ratio 3 between a and a p , the planetesimal 
and the perturber’s semi-major axis, respectively, and to average 
after over the two orbits. This permits us to perform an anal- 
ysis without any restriction on the eccentricities. The resulting 
Hamiltonian is given by Krymolowski & Mazeh (1999), Ford 
et al. (2000), or Beust & Dutrey (2006). 

To the lowest order in a (second order, quadrupolar), it yields 
a forced eccentricity ef. 


5 ae p 


Cf - t 


4 1—, 


2 ' 

P 


(2) 


This expression is given by Augereau & Papaloizou (2004) and 
Mustill & Wyatt (Eq. (8) of 2009). Note that this approach is 
only valid for low enough values of a to ensure a fast con- 
vergence of the expansion, that is, for orbits with significantly 
different sizes. It is also only valid far from mean-motion reso- 
nances. However, these resonances’ spatial extension along the 
semi-major axis (of the order of ~0. 1 AU) is typically two orders 
of magnitude smaller than the extent of the observed structures 
(of the order of ~ 10 AU), although when the particles are on ec- 
centric orbits, these resonances may span much wider ranges in 
terms of radial distance to the star than their span in semi-major 
axis may have hinted at. But in any case, the amount of mate- 
rial trapped in resonance can reasonably be assumed to be much 


3 u is such that a < 1, always, and thus a = c^/a if a p < a, and 
inversely, a = a/a p if a p > a. 
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smaller than the amount of non-resonant material, all the more 
so because we did not assume here that the planet has migrated, 
and thus excluded resonant capture during migration. Therefore, 
these structures are assumed to result from non-resonant mate- 
rial, and our approach is appropriate. Moreover, resonances were 
be treated in our A-body integrations, and were confirmed to be 
only important for limited parameter combinations (see Sect. 4). 

To derive the precession rate in the spherical harmonic ex- 
pansion case, we followed the method given by Mardling & Lin 
(2002). The variation rate for the Runge-Lenz vector of the orbit 
was computed, expanded to any given order, and integrated over 
one orbital period. 

After the numerical tests, there were fewer than two orders of 
magnitude between terms of the second and the fourth order (the 
third-order terms cancel out), therefore we retained terms up to 
fourth order in the spherical harmonic expansion, and averaged 
the resulting precession rate over the longitude of periastron. The 
precession rate is 

dm 3 n m p 3 Vl - e 2 

~dF = (1 - e 2 ) 372 

45 n m p . (4 + 3e 2 )(2 + 3e 2 p ) 

+ 256 M* a (1-e 2 ) 7 ' 2 

We then followed Wyatt et al. (1999) and assumed that the pre- 
cession timescale, r pre c, corresponds to the lower limit of the typ- 
ical dynamical timescale for setting a dynamical steady state: 

fprec = (dwd 7)Z’ (4) 

where a c is the typical semi-major axis of a particle in the ring. 

With this we analytically predict the effect of a perturber on a 
debris disk, that is, for a given set of values of e p and of the planet 
periastron q p , one can derive the precession rate (dtzr/d f) 0c corre- 
sponding to planetesimals that orbit at this distance. Conversely, 
one can set this dynamical timescale and the forced eccentric- 
ity for a particle with semi-major axis a c to correspond to those 
derived from observations of an eccentric debris disk, and thus 
initially estimate the perturber’s characteristics (see Fig. 2). 

However, the problem is more complex for real disks, which 
have a finite spatial extension, since these estimates depend on 
radial locations. To first order, it can be seen from Eqs. (2), to (4) 
that the forced eccentricity and the secular timescale scale as 


ei oc a, (5) 

and 

fprec ^ 7 ' (6) 

m p a 3 

We now define and a ou t as the inner and outer limits of the 
disk in the semi-major axis and define e^ m i n and ef >max as the 
minimum and maximum eccentricities induced across the disk 4 * . 
The minimum and maximum precession timescales, f pre c,mm 
and tprec.maxs are defined in the same manner. Then, using Eqs. (5) 
and (6), one obtains 

^f,max _ ^out 
£f,min ^in 

4 £f,min = £f,in and ef >max = £f, 0 ut for an inner perturber and vice versa 

for an outer one. 



160 200 260 300 360 

Perturber pericentrc ( AL') 


Fig. 2. Example colour map of the maximum induced eccentricity 2ef 
imposed by a planetary perturber on a particle with semi-major 
axis 100 AU and eccentricity e = 0 as a function of its periastron and 
eccentricity, as estimated from Eq. (2). The black line corresponds to 
a 2e f = 0.3 condition, which was set to mimic the condition for the disk 
of f 2 Reticuli. Note that it does not depend on the mass of the planet. 
The white lines show the parameters for which the typical timescale to 
reach a steady state at 100 AU is f prec = 1 Gyr, using Eq. (4). This 
timescale depends on the mass: m p = 0.1 Mj up (solid line), 1 Mj up 
(dashed line) and 2 Mj up (dotted line). For example, a perturber of 
mass 0.1 Mj up , periastron q v = 150 AU and eccentricity e p = 0.4 is 
expected to produce a significantly eccentric ring in shorter than 1 Gyr, 
although spiral patterns may remain since it can take several precession 
timescales for them to vanish, as was shown by Wyatt (2005). 


and 


fprec.max 
fprec. min 



( 8 ) 


It is easy to see from these equations that the secular preces- 
sion timescale spans a wide range of values across the disk. 
This means that making analytical predictions by setting the de- 
sired values for the forced eccentricity and the secular precession 
timescale for a particle with semi-major axis at the centre of the 
distribution suffers limitations when applied to an extended disk, 
especially concerning the timescale. 

Equations (7) and (8) can be rewritten using A a, the half 
width of the disk extent, along with ef, c and r prec ,c, the forced 
eccentricity and secular precession timescale at a c , respectively: 


e f 

,max/min 


a c ± Aa\ 

e f,c> 

a c I 


(9) 


and 

I a c ± Aa j 3 ^ 

fprec, max/min — I I fprec.c" (Id) 

\ a c ) 

As an example, we set a c = 100 AU, A a = 25 AU, 2e c = 0.3 
and fp re c,c = 1 Gyr. These values are close to those derived for 
the disk of £ 2 Reticuli (e > 0.3 and extent 70-120 AU: Eiroa 
et al. 2010, and Sect. 3 of the present work). One obtains 

/ 2cf im i n / max = 0.225—0.375 . .. 

\ fprec, min/max = 0.42-1.95 Gyr. 

In these conditions, the extent of the disk is not expected to affect 
the global eccentricity of the disk too much, that is, we expect 
to recover on average a global eccentricity corresponding to the 
forced eccentricity at a c after the steady state is reached. But 
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Fig. 3. Herschel / PACS images of f 2 Reticuli at 70, 100 and 160 microns from left to right. North is up and east is to the left. The insets in the 
bottom-left corners show the PSF. The north-west lobe is noted N-W, while the south-east lobe is noted S-E. 


the problem is that the extent of the disk strongly affects the 
timescale to reach this steady state. This is a limitation of the 
analytic approach that can be overcome by the use of numerical 
simulations. 

3. Numerical investigation: a typical set-up, 
the highly eccentric, old disk of g 2 Reticuli 

To move beyond the simplified analytical approach and explore 
the high-eccentricity case on Gyr timescales, we resorted to 
numerical tools. We placed ourselves in the frame of the re- 
stricted three-body problem, thatb is, one central star, a planet, 
and a mass-less planetesimal. The symplectic IV-body code 
SWIFT-RMVS of Levison & Duncan (1994) was used to inte- 
grate the evolution of a ring of planetesimals around a solar-mass 
star, over 1 Gyr. We used a typical timestep of ~ 1/20 of the short- 
est orbital period and ensured a conservation of energy with a 
typical error of ~ 10 -6 on relative energy. This code has a crucial 
advantage over an analytical approach: it is able to handle close 
encounters and scattering processes, along with the short-term 
variations of the planetesimals orbital elements, whereas these 
effects were ignored in the analytical approach, for which orbits 
were averaged, short-term variations were los,t and the approach 
is valid only for a <K 1, that is, far from close encounters. As we 
show in Sect. 4, the scattering events play a crucial role in the 
system’s evolution. 

For the sake of clarity, the f 2 Reticuli system was consid- 
ered as a proxy for a typical mature and significantly eccentric 
debris ring. We explore different planet-disk configurations, and 
produce synthetic images for comparison with Herschel/PACS 
observations. We show that the hypothesis of an eccentric debris 
disk around f 2 Reticuli is fully consistent with the observations. 

3.1. The Reticuli system 

The G1V solar-type star f 2 Reticuli (HR 1010, HIP 15371) 
(Eiroa et al. 2013) is located at 12 pc (van Leeuwen 2007), has 
a luminosity L* = 0.97 Lq, loggr = 4.51, and is ~2-3 Gyr 
old (Eiroa et al. 2013). It has a binary companion ^'Reticuli, a 
G2-4V (Gray et al. 2006; Torres et al. 2006) star located at a pro- 
jected distance of 37 1 3 AU from f 2 Reticuli (Mason et al. 200 1 ). 
Bayesian analysis by Shaya & Oiling (201 1) of the proper mo- 
tions of these stars indicates a very high (near 100%) probability 
that the pair are physically connected. 

The presence of dust around f 2 Reticuli has been probed 
with Spitzer (Trilling et al. 2008); the results suggest a ~ 150 K 
emission at ~4.3 AU. However, the angular resolution of 


Spitzer is limited, and the dust spatial distribution remained un- 
constrained. New observations with HerschelfPACS completed 
the spectral energy distribution (SED), providing the sugges- 
tion of an optically thin, ~40 K, emission at ~100 AU, with 
fractional luminosity T-dust/I-* ~ 10 -5 (Eiroa et al. 2010). 
Moreover, Herschel/PACS provided spatially resolved images 
of the dust thermal emission surrounding £ 2 Reticuli at 70 pm 
and 100 pm (Eiroa et al. 2010). We present here newly re- 
duced Herschel / PACS images (see Fig. 3). The images show a 
double-lobe feature, asymmetric both in position and brightness. 
Note that at 70 pm, the probability for alignment with a back- 
ground source within 10" is extremely low, namely 10 -3 (Eiroa 
et al. 2010). The disk is not resolved at Herschel/SPIRE wave- 
lengths: newly reduced images and star-disk flux measurements 
are presented in Appendix A. 

As suggested by Eiroa et al. (2010), the asymmetry revealed 
by Herschel / PACS in the disk of f 2 Reticuli can be interpreted as 
a ring-like elliptical structure with e > 0.3 seen close to edge-on 
and extending from ~70 to ~120 AU, which is fully consistent 
with the information derived from the SED (Eiroa et al. 2010). 
Alternatively, it might also be interpreted as two clumps from 
an over-density of dust and planetesimals. In Appendix B, we 
investigate the system inclination on the line of sight, a crucial 
parameter required for correctly interpreting the observed struc- 
tures. More precisely, we determined the star inclination and as- 
sumed that the disk and the star are coplanar. The 50% probabil- 
ity value is i = 65.5°, meaning that the system is highly inclined 
to the line of sight, which tends to support the eccentric-ring sce- 
nario. 

Without a doubt, this asymmetric structure provides evi- 
dence that “something” is dynamically sculpting the disk. This 
may be the stellar companion Reticuli or a (as yet undetected) 
planet. The latter hypothesis is fully compatible with radial ve- 
locity measurements of f 2 Reticuli, which suggest there is no 
Jupiter-mass (or larger) planet interior to ~5-10 AU (Mayor 
et al. 2003), but which put no contraints on a small planet or 
a Jupiter-like planet at larger radii. It is also compatible with 
growing observational evidence for planets at large orbital sep- 
aration, that is, several tens to a few hundreds of AU from their 
host star (see, e.g., Luhman et al. 2007; Kalas et al. 2008; Marois 
et al. 2008, 2010). 

Constraints from direct imaging are presented in Fig. 4 using 
the two evolutionary models COND 2003 (Baraffe et al. 2003) 
or BT-settl 2011 (Allard et al. 2011). Details of the reduction 
procedure are presented in Appendix C. These constraints were 
obtained from VLT/NaCo archival data taken in November 2010 
in the Ks band. These data do not provide constraints on 
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Distance to the star (arcsec) 
1.0 1.5 2.0 



Table 1 . Initial parameters for the planetesimals used in our (V-body 
simulations. 


Parameters 


Values 


Number of particles 
Semi-major axis (AU) 
Eccentricity 
Inclination (°) 


150000 

Amin — 20 , Umax — 140 
t?min 9 , t? m ax — 0.05 
i m i n — 3 , (max — 3 


Fig. 4. Detection limits set by direct imaging on the presence of brown 
dwarf/close binary between 1" and 2.5" in projected separation. 


companions beyond a projected distance of ~30 AU. The pres- 
ence of a brown dwarf within ~20 AU is still compatible with 
observations. 


the approach used by Chiang et al. 2009). The chaotic zone is 
defined as the region where mean-motion resonances overlap. 
The width of this zone depends on the mass ratio between the 
central star and the perturber, 

Aa _ fledge “ Aplanet _ ^ ^ 2/7 ^2) 


^planet 


^planet 


where /i = m p i anet /m* (Wisdom 1980; Duncan et al. 1989). 
Consequently, one can deduce the semi-major axis of a planet 
of a given mass that generates a disk inner edge at 70 AU: 


Aplanet — 


fledge 


1 + 1.5/j 2 ^ 


(13) 


3.2. Numerical set-up 

We considered a ring of 150000 mass-less planetesimals uni- 
formly distributed between 70 and 140 AU (except when 
specified otherwise) around a solar-mass host star, with initial 
eccentricities randomly distributed between 0 and 0.05, and ini- 
tial inclinations between ±3°, while the remaining angles (lon- 
gitudes of nodes and periastra) were randomly distributed be- 
tween 0 and 2n. These values are summarised in Table 1 . This 
reasonably well mimics the low eccentricities and inclinations 
expected at the end of the protoplanetary phase. The radial extent 
of the model disk was configured to closely match the observed 
properties of the disk around c 2 Reticuli. 

Using mass-less test particles removes self-gravity in the 
disk, as well as a back-reaction of the disk on the planet. In gen- 
eral, both of these phenomena are significant when the planet 
mass is similar to the disk mass. We have no mass estimate 
for the debris disk of f 2 Reticuli. However, a well studied case 
is the debris disk of Fomalhaut, whose mass was estimated to 
be ~3-20 M e (Wyatt & Dent 2002; Chiang et al. 2009). Since a 
debris disk loses material over time due to the combined effects 
of collisional evolution, Poynting-Robertson drag, and radiation 
pressure, and since Reticuli is much older than Fomalhaut 
(440 Myr; Mamajek 2012), it is reasonable to assume that the 
debris disk surrounding f 2 Reticuli does not contain more than 
a few Earth masses. In this case, it is also reasonable to assume 
that the disk self-gravity and back-reaction are negligible, and 
the planet will still be able to imprint on the disk structure if its 
mass is at least 0. 1 M Jup ~ 32 M s . 

We considered two planet-disk configurations, inside and 
outside the planetesimal belt, and performing parametrically 
explored their influent orbital elements. 

3.2.1. Inner perturber 

The first case is that of a planet interior to the initial ring. In this 
case, we considered that the inner edge of the disk is located at 
70 AU, and made the classical assumption that it is truncated 
by the chaotic zone generated by a coplanar planet (see, e.g.. 


We chose to perform a parametric exploration of the mass and 
eccentricity of the perturber, fixing its semi-major axis to the 
value deduced from the formula above. The values explored 
were 0.1 M Sup , 0.5 M Jup , and 1 M Sup for the mass and 0.2, 0.4 
and 0.6 for the eccentricity (see Table 2). Here the disk initial in- 
ner edge was fixed halfway between 70 AU and the planet semi- 
major axis. 

In addition, we also considered a more massive perturber 
such as a brown dwarf, located farther inside the system. 
This was motivated by the observation of a suspicious point 
source on observations in As-band data from the ESO archive 
(ID 086.C-0732(A); PI: Lohne, 71574), which disappeared after 
re-reduction of the data, and consequently, since its existence is 
controversial, we chose not to show these images. However, if it 
were confirmed, this point source was showing a 42 _M Jup brown 
dwarf located at a projected distance of 17.5 AU from the star. 
Therefore, we investigated the possibility for the presence of a 
brown dwarf in the inner parts of the system. For an inner edge 
to be produced at 70 AU, the planet’s semi-major axis need to 
be 43.8 AU. Since this constraint is less relevant for very ec- 
centric inner perturbers (see Sect. 4), we chose here to set this 
value to the perturber’s apastron and not its semi-major axis. We 
fixed its periastron to 17.5 AU, which led to an orbital eccentric- 
ity of 0.43 and a semi-major axis of 30.6 AU. With this orbit, 
analytical predictions indicate that the perturberexcites planetes- 
imals eccentricities up to 0.4. Therefore, this orbit was chosen 
for numerical tests (see Table 2). 

3.2.2. Outer perturber 

The second case considered is that of a planet exterior to the 
ring. There is indeed growing evidence for planets at large or- 
bital separations, that is, several tens to a few hundreds of AU 
away from their host star (see e.g. Luhman et al. 2007; Kalas 
et al. 2008; Marois et al. 2008, 2010), and the mass constraints 
set by direct imaging are loose given the age of the system (see 
Sect. 4.1.3). Therefore, we also investigated the ability of an ex- 
ternal perturber to shape a disk. 

We considered coplanar outer planetary companions and ex- 
plored the impact of the eccentricity, mass, and periastron on 
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Table 2. Summary of numerical experiments with an inner and outer 
perturber, and for a brown dwarf and the stellar binary compan- 
ion f 1 Reticuli. 


Inner perturbers 

ra p 

(M Jup ) fl p (AU) 

e v = 0.2 

e p = 0.4 

01 

II 

O 

05 

0.1 

63.2 

I 

III, A 

III 

0.5 

59.8 

I 

II, B 

IV, c 

1 

57.9 

I 

II 

IV 

Outer perturbers 

ra p 

(M, u p ) c h (AU) 

e p = 0.2 

Ol 

II 

o 

4^ 

Ol 

II 

p 

Os 


150 

I 

II, D 

V 

0.1 

200 

I 

V 

v, E 


250 

V 

V 

V 


150 

I 

II 

III, F 

1 

200 

1 

II 

II 


250 

I 

I 

II 


150 

I 

III, G 

III 

2 

200 

I 

I 

I 


250 

I 

I 

II 


Other perturbers 

Perturber 

m p Orbital parameters 

Result 

Brown dwarf 

42 M Jup a p = 30.6 AU; e p = 0.43 

II 

f 1 Reticuli 

IMo Op = 2046 AU; e v = 0.815 

I 


Notes. Description of the disk at 1 Gyr: I) steady state, e < 0.2; 
II) steady state, e > 0.2; III) scattered disk; IV) resonant pattern; V) spi- 
ral pattern. The example cases highlighted in more detail in Sect. 4 are 
labelled from A to G. 


the disk asymmetry. While an inner edge is in general consid- 
ered as a clue for the presence of inner perturbers, it is obviously 
more of a conjecture to assume that an outer edge is formed in 
the same manner, since a disk instrinsically has an outer limit. 
Therefore, the outer edge was not assumed here to be formed be- 
cause of resonance overlap, and the planet periastron was fixed 
instead, to ensure that it does not cross the disk. To explore a 
great variety of situations despite the CPU-time limitations, we 
considered a rough parameter space consisting of all the possi- 
ble combinations between masses m p = 0.1-l-2Mj up , perias- 
tra q p = 150-200-250 AU, and eccentricities e p = 0.2-0. 4-0. 6 
(see Table 2). 

Note that the perturbers were set on an initially eccentric 
orbit, which requires an explanation, because we assumed that 
the disk is initally symmetric. Indeed, this pictures a situation 
where the process that sets the perturber on its eccentric orbit 
leaves the disk unperturbed. This may seem rather unrealistic, 
even though some scenarios can be envisaged. For instance, an 
inner perturber might acquire its eccentricity via a planet-planet 
scattering event. However, this event should be such that a single 
perturber remains in the system. Otherwise, additional perturba- 
tions of a second planet would generate an orbital precession of 
the eccentric perturber, which in turn would not be able to sculpt 
the disk into an eccentric shape. This is compatible with obser- 
vational constraints, since as previously mentioned radial veloc- 
ity measurements rule out a massive perturber in the inner sys- 
tem. For an outer perturber, an eccentric outer binary companion 
may be able to generate these initial conditions (see Sect. 6). In 
summary, retrieving realistic initial conditions relies on a com- 
plete study of the perturbations induced on the disk for multiple 
scenarios and, most probably, an extensive parameter space ex- 
ploration. This is beyond the scope of the present paper, and will 


be the subject of future work; this motivated us to choose simple 
initial conditions. 


4. Numerical investigation: results 

We present here results obtained for the two disk-planet set-ups 
we considered: an inner and outer planet, and when the per- 
turber is the stellar companion l x Ret. For each case, we aimed 
to find the result that gives an eccentric disk compatible with 
observational constraints. 

Eiroa et al. (2010) derived a lower limit for the eccentricity of 
the disk in £ 2 Reticuli of 0.3. Therefore, given the uncertainties 
in the estimation of the disk global eccentricity we computed 
from our simulations, a disk global eccentricity lower than 0.2 
was discarded in our analysis. This global eccentricity was eval- 
uated considering the geometry of an ellipse: the offset 6 of the 
centre of symmetry from one of the foci of an ellipse is simply 
the product of its semi-major axis a by its eccentricity e, that 
is, 6 = ae. For a disk from our simulations, 6 can be obtained 
by calculating the centre of symmetry of the disk, using the po- 
sitions of the test particles in the heliocentric frame: 6 is the dis- 
tance of this centre of symmetry to the star. The disk semi-major 
axis a was determined as follows: the disk was divided into su- 
perimposed angular sectors of 3°. For each of these sectors, the 
radial distribution of the particles was fitted to a Gaussian. This 
provided the radial position of the maximum density for each 
angular sector, and thus a set of points defining the global shape 
of the disk. It is then straightforward to retrieve a from this set of 
points by seeking for the major axis, that is, the maximum dis- 
tance between two opposite points. Finally, the disk eccentricity 
is simply e = 6/a. 

4.1. Inner perturber 

We chose four illustrative results (see Table 2). For each, we 
show pole-on projections of the system at 1 Gyr in Fig. 7. We 
also summarise the results of our simulations in Table 2. 

The best candidates needed a significant orbital eccentric- 
ity of ~0.4. The example of Case B is shown in Fig. 7 (upper- 
right panel). However, scattering processes may be important, 
and studying them allowed us to place constraints on the mass 
of the perturber. 

4.1.1. Scattered disks 

Inner perturbers may lead to very significant scattering pro- 
cesses. Namely, they can fill the inner parts of the disk instead 
of producing a well-defined ring. Such effects are presented with 
Case A in Fig. 7 (upper-left panel). 

This effect appears for low-mass perturbers. These per- 
turbers do not scatter material efficiently enough. This material 
then tends to populate the inner parts of the system. 

As a consequence, there is a lower mass limit for inner per- 
turbers, and in the specific case of the f 2 Reticuli system, this 
lower limit is between 0.1 and 0.5 A/ Jup . However, one cannot 
exclude that another more massive planet produces scattering of 
the material, blowing it out and leaving an inner hole, therefore 
a more correct way to express this constraint would be that per- 
turbers with masses as low as 0.1 Mj up should be accompanied 
by another more massive planet to create this pattern. But this 
scenario presents difficulties: while this second planet must be 
massive enough to clear the inner parts of the system of its ma- 
terial, it also must have a limited dynamical effect on the first 
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Fig. 5. Semi-major axis versus eccentricity diagram of the disk at 1 Gyr 
for the Case C perturber. Overdensities of planetesimals at ~79 
and ~95 AU correspond to planetesimals in 3:2 and 2:1 mean-motion 
resonance with the perturber. 



planet that sculpts the disk: this second planet must be distant 
and not too massive for the orbit of the first planet to remain 
unperturbed. Otherwise, this orbit would precess and no longer 
lead to an eccentric pattern. It is not our purpose here to investi- 
gate this scenario, but based on the previous arguments, it would 
most probably work in a very limited parameter space. 


4.1.2. Resonant patterns 

Remarkably, for very eccentric (e p = 0.6) inner perturbers of 
mass between 0.5 and 1 M Jup , resonant clumpy structures may 
arise. An example is Case C, shown in Fig. 7 (bottom-left panel). 
In Fig. 5, we show a semi-major axis vs eccentricity diagram of 
the disk: it reveals two populations in resonance with the planet, 
namely the 3:2 and 2:1 mean-motion resonances, at ~79 AU 
and ~95 AU, respectively. These resonant populations are shown 
in Fig. 6. Typically, a population in 3:2 resonance is associ- 
ated with capture during outward planetary migration, as is the 
case in our own solar system (Malhotra 1993, 1995). But here, 
these structures is most probably appear because we used the 
chaotic-zone formula (Eq. (13)) to determine the perturber semi- 
major axis, which was derived for perturbers on circular orbits. 
Therefore, it is expected that this formula works less efficiently 
with increasing orbital eccentricity of the perturber: the result 
is that the planet digs into the disk, and consequently, planetesi- 
mals unprotected against close encounters by mean-motion reso- 
nance are scattered out, leaving the resonant structures apparent. 
This is supported by the fact that these resonant structures dis- 
appear when the constraint given by the chaotic-zone formula 
is applied to the perturber’s apastron instead of its semi-major 
axis, as for an inner brown-dwarf-type companion. Interestingly, 
the observation of these resonant structures in a system may pro- 
vide other clues on the dynamical history of a perturber than an 
outward migration: it may mean that the planet was originally 
shaping the inner edge of the disk before it was set on an eccen- 
tric orbit. However, these are thin structures, and when the disk 
is seen close to edge-on, as is the case for c} Reticuli, they would 
most probably be hidden by the non-resonant population. 


4.1.3. Brown dwarf 

Additionally, we investigated the possibility for the presence of 
a brown dwarf in the inner parts of the system, on an orbit such 
that the perturber excites planetesimal eccentricities up to 0.4. 
The disk at 1 Gyr is shown in Fig. 7 (bottom-right panel). Its 
global eccentricity is ~0.2-0.25, which shows that very massive 
perturbers in the inner parts of the system can create the desired 
pattern. 

4.2. Outer perturber 

We chose four illustrative results (see Table 2) and show pole-on 
projections of the system at 1 Gyr in Fig. 8 for each, along with 
the evolution of the disk offset on 1 Gyr in Fig. 9. The results of 
our simulations are summarised in Table 2. 

4.2.1. Good candidates 

Case D is shown in Fig. 8 (upper left panel). The corresponding 
evolution of the offset clearly shows that the disk is at steady 
state and significantly eccentric (see Fig. 9). The best candi- 
dates must have significant eccentricities: none of our 0.2 ec- 
centric perturbers creates the desired eccentric structure, even 
in the limit case predicted analytically, where a 0.1 M Jup per- 
turber has a pericentre q p = 150 AU and eccentricity e p = 0.2 
(see Fig. 2). As a consequence, the best outer candidates have 
eccentricities ~0.4-0.6. 

4.2.2. Spiral patterns 

Case E is shown in Fig. 8 (upper-right panel). The system shows 
a spiral pattern at 1 Gyr, which, according to analytical pre- 
dictions, corresponds to one precession timescale (see Fig. 2). 
Thus, our A-body integrations confirm what was noted by Wyatt 
(2005): the analytical formula appears to be a lower limit and 
several precession timescales are sometimes necessary for spiral 
patterns to vanish. The effect of spirals can also be seen in the 
evolution of the disk offset, which oscillates (see Fig. 9). 

If we had observational proof that the disk of Reticuli has 
reached steady state and contains no spiral pattern, the results of 
our numerical investigation would allow us to place a lower mass 
limit of 0. 1 M Jnp on an outer perturber in a range of periastra 
from 150 to 250 AU, based on a dynamical timescales criterion, 
otherwise it takes longer than 1 Gyr to generate a steady-state 
eccentric disk. This limit still holds for larger periastra than the 
range explored, since dynamical timescales increase with dis- 
tance. But because of the slightly edge-on orientation of the disk 
and the resolution of the Herschel / PACS images, it is extremely 
difficult to discard a possible spiral pattern in this disk, and no 
lower mass limit can be clearly defined for an outer perturber. 

4.2.3. Scattered disks 

Outer perturbers may lead to very significant scattering pro- 
cesses. We present such effects in Fig. 8 (lower panels), where 
Case F and Case G are considered. 

These processes are even more significant when the mass 
of the perturber increases. Indeed, as physically expected, more 
massive perturbers tend to scatter small bodies more efficiently. 
The distance to the disk plays a major part in this effect too, since 
close perturbers also tend to scatter more material. Additionally, 
when a perturber is on an eccentric orbit, it approaches even 
closer to the disk. In the most dramatic cases, the disk is 
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Fig. 6. Views from above the plane of disks at 1 Gyr in Case C, i.e., where resonant patterns appear. Planetesimals in 3:2 (left) and 2:1 (right) 
mean-motion resonance with the perturber are shown in red. 
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Fig. 7. Example views from above the plane of disks at 1 Gyr, under the influence of inner perturbers. Top left', the perturber scatters material in the 
inner parts of the system (Case A). Right: one of the best candidates (Case B). Bottom left: resonant patterns appear (Case C). Right: a brown-dwarf 
perturber can be a good candidate. 
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Fig. 8. Example views from above the disk plane at 1 Gyr, under the influence of outer perturbers. Top left : one of the best candidates (Case D). 
Right: spirals are still apparent (Case E). Bottom left: one might question the contribution of the scattered material to the disk emission (Case F). 
Right: scattering processes tend to destroy the structure (Case G). 


completely destroyed. Consequently, for a given distance to the 
disk, there is an upper limit to the mass of an outer companion. 

For Case F, one might question the contribution of the scat- 
tered inner material to the emission of the disk, and whether it 
would be visible on resolved images. In this case, the potential 
visibility of material on real observations relies on the sensitivity 
and resolution of the instrument used for these observations, as 
well as on the distance of the object, the radiative properties of 
the material itself and the quantity of light it receives, that is, on 
the host star properties. Therefore, only the production of syn- 
thetic images for direct comparison with observations can reveal 
whether this material is apparent or not, and refine constraints 
on the potential perturbers at work in the system. Additionally, 
the evolution of the offset clearly suggests that the asymmetry 
relaxes asymptotically to a low but non-zero value, and indeed, 
the apparent ring structure shows little eccentricity. 

More specifically, our results allow us to place an upper mass 
limit of 2 A/j U p on an outer perturber in a range of periastra 
from 150 to 250 AU for the f 1 Reticuli system. However, this 


upper mass limit is expected to increase for periastra higher than 
those explored, since scattering processes are expected to be less 
efficient for a given mass when the companion is farther away 
from the disk. 


4.3. Stellar binary companion 

We investigated the influence of the binary companion 
if 1 Reticuli on the debris disk surrounding f 2 Reticuli. The aim 
was to determine whether it alone could generate the observed 
asymmetry. If this were to be the case, the companion must be 
on an eccentric orbit. The only observational constraint available 
so far on the relative orbit of the binary system is the projected 
distance between the two stars, namely 3713 AU (Mason et al. 
2001). Consequently, we cannot exclude the possibility that its 
orbit is eccentric, with a present-day location near apastron. We 
investigated whether a binary companion on an eccentric orbit 
at this distance might account alone for the observed structure 
of the disk, which would mean an elliptic ring with a minimum 
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Fig. 9. Example of the time evolution on 1 Gyr of the disk offset coordinates ( X (solid line) and Y (dashed line)) for outer perturbers. Top left : 
one of the best candidates (Case D). Because of the orientation of the perturber orbit, where its semi-major axis is along the X-axis with positive 

periastron, this corresponds to a negative X-offset AX 30 AU and a zero F-offset for a disk located at a c ~ 100 AU. Top right: the oscillations 

of the offset coordinates reveal the spiral-winding regime (Case E). Bottom: the offset seems to relax; left: case F; right: case G. 
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Fig. 10. Semi-major axis vs. eccentricity diagram (left) and pole-on projection (right) of a disk of planetesimals after 2 Gyr under the influence of 
a one-solar-mass star of semi-major axis 2046 AU on an orbit of eccentricity 0.815. The disk is truncated at ~80 AU and nearly circular. It has no 
cleai' offset: its centre of symmetry is offset by ~6 AU from the star along the perturber major axis. With a disk mean semi-major axis of ~70 AU, 
this gives a global eccentricity of ~0.08. 


global eccentricity of 0.3, without any constraint on the binary 
eccentricity. We used analytical and numerical methods. 

We first examined the forced secular eccentricity applied to 
the planetesimals by an eccentric binary companion, which was 
assumed to be coplanar. The debris disk surrounding f 2 Reticuli 
is approximately centred on a = 100 AU, and the binary per- 
turber is at 3713 AU from f 2 Reticuli. This is only a projected 
distance, not a semi-major axis. We assumed that the binary 


companion is currently located at a distance r* from f 2 Reticuli. 
The equation of its orbit around f 1 Reticuli reads 


'★ — ; v, A ■ 7 

1 + e* cos u* 

where u*, a* and e* are the binary companion current true 
anomaly, semi-major axis, and orbital eccentricity, respectively. 
The observed distance d = 3713 AU is related to r* 
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by d = r* cos if/, where if/ is a projection angle. This gives 
d( l+e*cosu*) 

a * — 2\ i (15) 

(1 - e+) cost// 

From this result Eqs. (15) and (2), where a p and e p are substituted 
by a* and e*, 

5 a e* cos i/r 

2e f ---— (16) 

2 a 1 + e* cos i>* 

It is clear from this equation that the highest possible ef values 
ware obtained for cosu* = -1 (binary currently at apastron) 
and cosi/r = 1 (no projection factor). With these assumptions, 
one derives 

2e f ,max - 0.068-—^— — (17) 

1 - e* 

For 2e f , max to reach at least 0.3, e* > 0.815 is required. 

This seems highly eccentric and very unlikely at first sight. 
Yet, Duquennoy & Mayor (1991, see their Fig. 6b) have 
shown that almost 25% of binaries with orbital periods longer 
than 10 3 days have orbital eccentrities e* = 0.825 ± 0.075. In the 
present case, d = 3713AUande* = 0.815 lead to a* = 2046 AU 
and an orbital period T* ~ 10 5 yr. Therefore it is possible 
that f 1 Reticuli is on an eccentric orbit, if not a highly eccentric 
one. However, the derived orbit should also have an apoastron 
value of q t = 379 AU and one might question the disk sur- 
vival at ~70-120 AU with a stellar- type perturber approaching 
so close to the system. 

Equation (1) of Holman & Wiegert (1999) gives the critical 
semi-major axis a cr it for orbital stability around a star perturbed 
by a binary. This is 

a crn = [(0.464 ± 0.006) + (-0.380 ±0.010)/z 

+(-0.631 ± 0.034)e* + (0.586 ± 0.061)yre* 

+(0.150 + 0.041)4 + (-0.198 + 0.074)^4]a*> (18) 

where = m*/(m^- 2 Reticuli + m*) is the star mass ratio of 
value 1/2 if we assume here wi* = mg. Reticuli = 1 Af®. Material 
with a > dcrit will be on an unstable orbit, and most probably 
scattered out of the system. 

Applying this to our case leads to a cr i t = 66/j?g 6 AU. 
Uncertainties on a cr j t are very large, and this result shows that 
within uncertainties, the disk could exist at the observed dis- 
tances, or, in contrast, not exist at all. Therefore, this orbit was 
tested numerically. We considered a ring of 150000 mass-less 
planetesimals uniformly distributed between 70 and 140 AU 
from their solar-mass host star, with proper initial eccentricities 
randomly distributed between 0 and 0.05, and initial inclinations 
between ±3°. 

The test particles are perturbed by another solar-mass star 
in orbit around the primary with a semi-major axis 2046 AU 
and eccentricity 0.815, coplanar to the disk. We used the 
SWIFT-RMVS /V-body symplectic code of Levison & Duncan 
(1994) to compute their orbital evolution for 2 Gyr. 

Our results are consistent with the predictions of Holman & 
Wiegert (1999): the disk is truncated at ~80 AU (see Fig. 10, 
left panel), which is incompatible with observational constraints 
of the debris disk of £ 2 Reticuli, since it radially extends up 
to ~ 120 AU. 

Moreover, at 2 Gyr, the disk has a highly symmetric shape 
and no clear offset from the star (see Fig. 10, right panel). These 
results suggest that the binary companion on an eccentric orbit 
can not account for the disk structure, and that an unseen eccen- 
tric companion is more likely responsible for shaping the disk. 


Another possible way for a stellar binary companion to gen- 
erate high-eccentricity orbiting bodies around a primary is the 
Kozai mechanism (Kozai 1962). This would occur if the disk 
and the binary orbit were mutually inclined by more than ~39°. 
In that case, the orbit of a particle in the disk would suffer 
coupled modulations in inclination and eccentricity, that is, a 
particle would periodically switch from a highly eccentric or- 
bit, coplanar with the binary companion’s orbit, to a circular 
orbit, highly inclined with respect to the orbit of the binary 
companion. However, we can discard this mechanism as a pos- 
sible explanation for the eccentric global structure of the disk 
of 4 Reticuli: indeed, the Kozai Hamiltonian is invariant by ro- 
tation, meaning that, if this mechanism is able to excite planetes- 
imal eccentricities to high values, their longitudes of periastra 
remain uniformely distributed, while they should be preferen- 
tially aligned for an eccentric and offset disk to be generated. 


5. Synthetic images 

In this section, we use the results of our best-fit A-body simula- 
tions to produce a realistic dust population and synthetic images 
for comparison with Herschel / PACS observations. The proce- 
dure followed to create synthetic images is straightforward: a 
population of dust was created from the position of the parent 
planetesimals. 

The main difference between dust particles and planetesi- 
mals is that the former are small enough to be affected by stellar 
radiation pressure. Radiation pressure is usually described for a 
particle by its constant ratio j3 to stellar gravity. The dust particles 
are assumed to be released by planetesimals, which feel no radi- 
ation pressure. Hence the daugther particles assume an orbit that 
is very different from that of their parent bodies. It is well known 
that if the parent bodies move on circular orbits, the dust parti- 
cles are unbound from the star as soon as (3 > 0.5. In our case, 
however, dust particles may be released by planetesimals orbit- 
ing on more or less eccentric orbits, which may slightly change 
this threshold. Because planetesimal eccentricities are expected 
to be moderate on average, y3 = 0.5 can nevertheless be consid- 
ered as a reasonable approximation. 

Small grains are released from seed planetesimal positions 
at the planetesimal velocity, and are then spread along the orbits 
determined by these initial conditions and their y3 value. We are 
aware that this simple procedure cannot accurately evaluate the 
spatial distribution of the smallest grains. To do this, complex 
models such as the DyCoSS code of Thebault (2012), the CGA 
of Stark & Kuchner (2009), or the LIDT-DD ode by Krai et al. 
(2013) have to be used to evaluate the complex interplay between 
the rate at which grains are collisionally produced from parent 
planetesimals, the time they spend (because of their highly ec- 
centric orbits) in empty collisionally inactive regions and the rate 
at which they can be affected or even ejected by close encounters 
with the perturbing planet (see Sect. 4 of Thebault et al. 2012), 
not to mention the Poynting-Robertson drag these small grains 
are subject to. 

However, this caveat is acceptable for the present problem, 
because the role played by small micron-sized grains close to 
the blow-out size is very minor at wavelengths >70 yim, so that 
our synthetic images are not strongly affected by errors in their 
spatial distribution. 

We set the dust grain sizes range from 0.5 yim to 1 mm, with a 
classical Dohnanyi (1969) power-law distribution (index -3.5), 
which covers the [5 distribution from 0 to 0.5 well, since this 
parameter depends on grain size. 
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Table 3. Stellar and disk fluxes at 70, 100, and 160 fim (Eiroa et al. 
2010 ). 


Wavelength (pm) 

Stellar flux (mjy) 

Disk flux (mjy) 

70 

24.9 ± 0.8 

8.9 ±0.8 

100 

13.4 ± 1.0 

13.5+1.0 

160 

19.4 ± 

1.5 


~4.7 

-14.7 


Notes. The fluxes at 70 and 100 gm, along with the total flux of the star- 
disk system at 160 pm are PACS measurements. The individual star and 
disk fluxes at 160 pm result from predictions. 

Table 4. Disk position angle observed with Herschel/PACS , disk po- 
sition angle in our synthetic images, and telescope orientation during 
Herschel/PACS observations and PSF. 


Wavelength (pm) PA(°) 


Observed disk 

- 

110 

Simulated disk 

- 

110 

PSF 

70/160 

127 

PSF 

100/160 

127 

i f Reticuli 

70/160 

281 

t ( 2 Reticuli 

100/160 

281 


Their emission was computed using the radiative transfer 
code GRaTeR (see, e.g., Lebreton et al. 2012). To do this, the 
following parameters are required: distance of the star (12 pc), 
magnitude in band V (V = 5.24), and total luminosity 0.96 Lq. 

Because the disk is optically thin, its mass is linearly linked 
to the flux emission intensity, and it can be easily scaled to fit the 
observed intensity (see Table 3). The mass needed for the disk to 
produce a flux as observed on Herschel/PACS will vary with the 
dust grain composition, and thus its density. But because we have 
no constraints on the dust composition, astrosilicate grains were 
used (Draine 2003), and the mass of the disk was simply scaled 
to obtain intensities compatible with observational constraints 
for a given wavelength. 

Thermal emission images were produced with a resolu- 
tion of T'/pixel at 70 and 100 /mi and 2"/pixel at 160 pm. 
Before convolving these images with the point spread function 
(PSF), the star was added at the central pixel with a flux inten- 
sity matching the predicted stellar photosphere flux density in 
each waveband. The position angles of the disk observed with 
Herschel/PACS and of the disk in our synthetic images, as well 
as the orientation of the telescope during the observations were 
taken into account (see Table 4). Our purpose was to match the 
observations. 

We chose the simulation that led to a clear and signif- 
icantly eccentric disk at 1 Gyr, namely our Case A (see 
Table 2, and upper-left panels of Figs. 8 and 9) to reproduce 
the Herschel/PACS image at 100 yum. The mass of the disk 
was scaled so that a total flux of 13.5 mjy as observed with 
Herschel/PACS at 100 pm was spread all over the disk. The 
last parameter needed is the system inclination. Our best fit 
gives 65.5° (see Appendix B). 

We present an unconvolved and convolved image of this 
disk in Fig. 1 1 . The disk offset is clearly visible in unconvolved 
images. However, there seem to be no difference between the 
symmetric and asymmetric state after the images are convolved 
because of the contrast between the star and the disk emission. 
Moreover, the disk flux per pixel is one order of magnitude lower 
than the fluxes observed in Herschel/PACS images. 


This is not surprising, however: to estimate the total disk 
flux, the flux was measured in a small region of the disk before 
applying aperture correction. Even if the correct aperture cor- 
rection for a point source were used, it would always be a lower 
limit, and the total disk flux would be underestimated. The par- 
ent ring may also be narrower, which would increase the flux per 
pixel. 

Therefore, we investigated the impact of the width and of 
the total flux of the disk on the features visible with PACS 
and produced convolved images of a dusty disk produced by 
an asymmetric eccentric parent ring of diverse total fluxes (1, 
2, and 5 times the flux measured by HerschelfPACS) and of di- 
verse widths (semi-major axis centred on 100 AU, widths 5, 10, 
and 20 AU). To do this, particles from a range of semi-major 
axes from our /V-body simulation output were selected. An in- 
clination angle of 65° was chosen, which is the most probable 
inclination derived for the system. 

The fluxes per pixel recovered in convolved images with a 
disk five times more massive than the mass initially derived from 
observations match the observations better. This provides a bet- 
ter constraint on the disk mass and total flux. We show this exam- 
ple in Fig. 12 (left panel): the asymmetric double-lobe structure 
now appears clearly. 

It is worth noting that the width of the parent ring has a more 
limited influence on the flux per pixel than the mass of the disk. 
But it has an influence on the appearance of the disk: a narrow 
parent ring (<10 AU) leads to a resolved apastron lobe, which is 
more consistent with Herschel/PACS images, although this lobe 
does not seem to be located as far from the star as it is in the 
Herschel/PACS images. 

Therefore, we also investigated the location of the disk by 
producing convolved images of a ring of dust produced by a nar- 
row eccentric parent ring of width 5 AU and semi-major axis 
centred on 120, 130, and 140 AU. The disk total flux was set to 
be five times the flux derived front observations. The inclination 
angle was here again set to 65°. As expected, the farther away 
the disk is located, the clearer the lobes appear (see Fig. 12). 

Finally, with all these insights, we conclude that the hypoth- 
esis of an eccentric dusty disk around Reticuli is indeed com- 
patible with Herschel/PACS images, provided that the dust is 
produced by a narrow parent ring with width less than 10 AU 
and located slightly farther away than derived by Eiroa et al. 
(2010), such that its a semi-major axis distribution is centred be- 
tween 120 and 140 AU. 

This slightly changes the constraints derived on potential 
perturbers, but the forced eccentricity only depends on the ratio 
between the planet and planetesimal semi-major axis, and in a 
linear way at lowest-order approximation. This means that in 
this approximation the constraints can be completely scaled in 
a linear way, that is, the potential semi-major axis for planets 
must also be increased by 20^10% and the disk needs to be 
centred on 120-140 AU, while the constraint on the perturber 
eccentricity (e > 0.3) remains identical. 

6. Discussion and conclusion 

We used <T 2 Reticuli as an example to discuss the shaping of a 
disk into an eccentric ring on Gyr timescales. 

We showed that eccentric patterns in debris disks can be 
maintained on Gyr timescales, but also that eccentric perturbers 
can produce patterns other than eccentric rings. 

The general results of our simulations show that both inner 
and outer perturbers can generate extremely significant scatter- 
ing processes. They can cause a disk to adopt structures that 
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Fig. 11. Synthetic images at A = 100 //m, unconvolved disk (no star) (left) and convolved disk and star image with the PSF (right). Disk at 1 Gyr 
in Case A (see Table 2), seen with inclination 65°. The star is at the centre of the image, and in the convolved image, the flux scale is set to match 
that in the Herschel / PACS image. 
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Fig. 12. Synthetic images of Case A (see Table 2) at 100 //m and convolved with the PSF. The parent ring is centred on 100 AU (left) and 130 AU 
(right), has a width of 5 AU and is seen with an inclination of 65°. The disk total flux is five times the flux measured by Herschel/ PACS. The star 
is at the centre of the image and the flux scale is set to match that of Herschel / PACS image. 


show no clear elliptic ring: the inner part of the disk is either 
filled or the structure is destroyed. These scattering processes en- 
danger the survival of an eccentric ring, and investigating these 
processes with numerical experiments allowed us to put con- 
straints on potential perturbers. From these constraints, we de- 
rived an upper-mass limit for outer perturbers in a certain range 
of periastra and a lower-mass limit for inner perturbers. 

Moreover, the timescale for spiral structures to vanish is 
longer with smaller-mass perturbers, and thus, investigating this 
timescale with numerical experiments permitted us to place a 
lower-mass limit on perturbers, provided that spiral structures in 
a disk can definitely be ruled out from observations. The role 
of numerical experiments is crucial here, since the analytical 
timescale at the disk centre of distribution underestimates the 
effective timescale in a real extended disk. 


The offset of the disk centre with respect to the star is mostly 
stable; however, we note that it seems to relax very slowly in 
rare cases. While the former evolution is characteristic for peri- 
centre glow dynamics, the latter one is surprising. The relaxation 
of the eccentric structure is not expected in the first-order secu- 
lar analysis described by Wyatt et al. (1999) and Wyatt (2005). It 
may be the result of higher-order terms that have been neglected 
in the analytical study, and is more probably caused by erratic 
short-term variations of the planetesimals’ semi-major axes due 
to moderately distant approaches to the planet. These effects, 
which can lead to scattering of the planetesimals, are elimi- 
nated in the analytical averaging process of the perturbations, 
and thus cannot be predicted analytically in the secular approxi- 
mations used here. A more detailed study of this relaxation phe- 
nomenon of pericentre glow structures is nevertheless beyond 
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the scope of the present paper and will be the purpose of future 
work. 

Transient spiral structures, filled inner holes, sparsely popu- 
lated scattered disks, and resonant clumpy structures are all pos- 
sible results when an eccentric perturber acts on a debris disk. 
They can be shown with numerical simulations, but more im- 
portantly, used to derive constraints on a perturber in a system 
(mass and eccentricity). Therefore, we provided a method for 
investigating and modelling eccentric ring structures based on a 
complementary analytic and numerical approach, where one can 
derive potential orbits from analytics and test them numerically 
using /V-body codes. This method can be easily applied to other 
systems and is expected to be useful in the near future. 

Indeed, Kaib et al. (2013) have pointed out that wide bi- 
nary star systems, that is, systems with separations greater than 
1000 AU, can produce eccentric planets around a primary star 
on Gyr timescales. This is due to Galactic tides and passing 
star perturbations, which are able, sooner or later, to set the sec- 
ondary star on a highly eccentric orbit. The proportion of wide 
binary systems is by no means negligible (~50%, Duquennoy & 
Mayor 1991), and although debris disks which are several Gyr 
old are faint and difficult to detect, this will be overcome with the 
unique capabilities of ALMA, JWST, and SPICA. Therefore, old 
eccentric patterns in debris disks are expected to be commonly 
observed in the future. 

The Reticuli disk is one such example of such a Gyr-old 
eccentric debris disk. Moreover, C Reticuli is part of a wide bi- 
nary star system, which may provide an explanation for the pres- 
ence of an eccentric perturber around l 2 Reticuli. We showed 
that the binary companion cannot be directly responsible for the 
eccentric ring structure, and also that the asymmetry is instead 
caused by a closer companion, either interior or exterior to the 
disk. In all cases, the eccentric companion is expected to have an 
eccentricity e > 0.3 to produce this pattern. 

Investigation of the disk structure generated by the scattering 
processes provided an upper-mass limit of 2 M Jup for an outer 
perturber located in a range of periastra 150-250 AU, whereas a 
lower-mass limit of 0.1 Mj up is associated with inner perturbers 
in the f 2 Reticuli system. 

By producing synthetic images, we showed that the origi- 
nal interpretation of the double-lobed feature around f 2 Reticuli, 
that is, the observed eccentric ring e > 0.3, is clearly supported, 
although the disk is located slightly farther away (20-40%) than 
originally derived. Moreover, we found that the dusty disk is 
probably created by a narrow parent ring (width <10 AU), which 
should have a slight inclination with the line of sight, compatible 
with the most probable inclination derived for the system, and it 
also has a significantly higher flux than that estimated from the 
HerschelfPACS measurements (at least five times). 
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Appendix A: SPIRE images of £ 2 Reticuli 



Fig. A.l. Herschel/SPIRE images of f 2 Reticuli at 250, 350, and 500 /rm ( left to right). The images were reduced using the standard reduction 
scripts on HIPE, version 8.2, and SPIRE CAL 8.1. Image orientation is north up, east to the left. The pixel scales are 6", 10", and 14" at 250, 350, 
and 500 /rm, respectively. The SPIRE beam FWHM in each band is denoted by the white circle in the bottom-left comer of each image. 


Table A.l. Predicted stellar fluxes and SPIRE flux measurements of the star-disk system at 250, 350, and 500 /rm (Eiroa et al. 2013). 


Wavelength (/rm) 

Stellar flux (mJy) 

Disk+star flux (mJy) 

250 

2.03 ± 0.03 

59.72 ± 6.70 

350 

1.04 ±0.02 

24.68 ± 6.89 

500 

0.51+0.10 

20.29 ± 7.66 
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Table B.l. Stellar properties of f 2 Reticuli. 
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Stellar property 

Value 

Reference 

(B-V) 

0.60 

1 

R* (Ro) 

-0.965 ± 0.05 

2 

U (L e ) 

0.97 

3 

T* (K) 

5851 

3 

u rot sin i (km s 1 ) 

2.7 ± 0.3 

4 


-4.79 ± 0.03 

5 

log(Ro) 

-0.185 + 0.085 

2 

r c (days) 

-9.10 

2 

P rot (days) 

-14.20 + 2.75 

2 

v mt (kms-1) 

-3.42 ± 0.66 

2 

in 

~65.5 +2 ^ 

2 


References. (1) Johnson et al. (1966); (2) this study; (3) Eiroa et al. 
(2013); (4) Reiners & Schmitt (2003); (5) Henry et al. (1996). 



Fig. B.l. Possible inclination angles (taken from the pole) 
for f 2 Reticuli as a function of u rot ranging between -2.76 
and 4.08 kms -1 , where i is computed using the values of 
e r otsin(i) = 2.7 ± 0.3 kms -1 (Reiners & Schmitt 2003). The acceptable 
zone for the inclination is grey-shaded. 


Observations of the debris disk surrounding £ 2 Reticuli re- 
veal a double-lobed asymmetric feature. The inclination of this 
system relative to the line of sight is a key parameter for correctly 
interpreting the observations. If the system is seen pole-on, one 
would expect the observed feature to be the signature of reso- 
nant clumps, whereas an eccentric ring signature would be more 
plausible if the system were observed edge-on. 

In general, observations suggest that stellar and disk inclina- 
tions are aligned (Watson et al. 2011; Guilloteau et al. 2011). 
Under this assumption, one can estimate the disk inclination 
from the observed stellar inclination. Consequently, we aim here 
to measure the star’s inclination i, that is, the angle formed by its 
rotation axis with respect to the line of sight. With this conven- 
tion, the system is seen pole-on if i = 0°, and edge-on if i = 90°. 

The method used requires knowing the following stellar 
properties: the colour index ( B - V), the radius R*, the pro- 
jected rotational velocity, v mt sin(i), and finally R[j K , an activity 
indicator defined as F' HK /crT 4 , where F^ K is the chromospheric 
flux in the H and K lines of Call, and T t is the effective star 
temperature. These properties for c 1 Reticuli are summarised in 
Table B.l. 

We first used the activity/rotation diagram built by Noyes 
et al. (1984), which plots log (R’ HK ) versus log (Ro) and shows 
a relationship between these two quantities for late-type stars. 
Ro = Rrot/rc is the Rossby number, P mt is the rotational pe- 
riod of the star, and r c a model-dependent typical convective 
time, called the turnover time. Using Fig. 6b of Noyes et al. 
(1984) and the observed value of log(RJj K ) = -4.79 found 
by Henry et al. (1996) for t 2 Reticuli allows us to estimate 
log(R ro t/Tc)~ 0.185 ±0.085. 

Then, using Eq. (4) of Noyes et al. (1984), where x is defined 
with the star colour index (B - V) by x = I — (B - V), one can 
estimate r c : 

, \ 1.362 - 0.166.V+ 0.025x 2 - 5.323x 3 ,x>0 . , 

l0gTc = { 1.362 - 0.14* ,*<0 ■ (B ' 1} 

f 2 Reticuli has a colour index (B - V) = 0.60 (Johnson et al. 
1966), which gives r c = 9.10 days. The corresponding range of 
possible rotation periods is R rot ~ 14.20 ± 2.75 days. 

Since the equatorial rotation velocity is defined as u rot = 
2 7rR*/R rot , knowing the stellar radius R* allows us to obtain a 
range of possible values for v TOt . Using T* and Lb 0 i and cor- 
rections prescribed by Masana et al. (2006), we found that for 


ip Reticuli, R* = 0.965 R 0 ± 0.05. The corresponding value of 
equatorial velocity is -3.42 ± 0.66 km s -1 . 

We compared this with the independent measurement of 
v rot sin i by Reiners & Schmitt (2003). They found v tot sin(i) = 
2.7 ± 0.3 kms -1 , which combined with the calculated rotational 
velocity, allows us to estimate the stellar inclination. We found 
that the inclination can range from 34° to 90° (see Fig. B.l). This 
is very consistent with an inclined disk. 

However, two angles are required to fully constrain the stel- 
lar rotation and disk axis. Therefore a degree of freedom remains 
and different orientations may lead to the same inclination i. 
Namely, the range of possible orientations leading to a same in- 
clination is the set of axes describing a solid angle 2n sin i about 
the line of sight. But the number of axes that lead to the same 
inclination increases with i, since it follows a sin i distribution. 
This means that the inclinations in a range [34°; 90°] are not 
equiprobable. The probability for a given inclination to be be- 
tween i and i + di, provided it is in the range [34°; 90°], can thus 
be written 


d P(i) 


sin idi 



(B.2) 


and the probability to find inclinations between i\ and i 2 in the 
range [34°; 90°] is 


P sin idi 

P(ie\i r ,i 2 \) = (B.3) 

J34 0 Sln/d ' 

Applying this to the case of f 2 Reticuli, among the possible 
ranges of [34°; 90°], we have an -50% chance that the observed 
inclination is in the range [65.5°; 90°]. Thus the system incli- 
nation is i = 65.5°3 24 5 5 . Consequently, the disk is more probably 
seen almost edge-on, with a pure edge-on configuration not ruled 
out. 

Because of the large uncertainties, this constraint does not 
really allow us to conclude with absolute confidence whether 
the disk exhibits an eccentric ring or resonant clumps. However, 
resonant structures are in general thin structures that tend to be 
hidden by non-resonant bodies and are difficult to detect, even in 
case of pole-on observations (Reche et al. 2008). This argument 
clearly supports the interpretation of Eiroa et al. (2010), who 
favoured an eccentric ring structure with e > 0.3 seen edge-on 
and extending from ~70 to ~ 120 AU. 
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Appendix C: Constraints on £ 2 Ret set by direct 
imaging 

VLT/NaCo (Lenzen et al. 2003; Rousset et al. 2003) Ifs-band 
data were retrieved from the ESO archive (ID 086.C-0732(A); 
PI: Lohne, 71574). Two epochs were available in August 2010 
and November 2010, the former missing photometric calibra- 
tion therefore only the latter was used to set detection limits on 
the presence of bound companions. Nevertheless, both data sets 
were reduced and no companion was detected. The data from 
November 11, 2010 were obtained in field stabilized-mode with 
five manual offsets of the derotator to simulate field rotation, 
with the S27 camera providing a pixel scale of 27 mas/pixel. 
Twenty image cubes with a DITxNDIT of 1.5 s x 42 were ob- 
tained for a total observing time on target of 2 1 min. The semi- 
transparent mask C_0.7_sep_10 with a diameter of 0.7" and a 
central transmission of 3.5 x 10 -3 was used. Each individual im- 
age was bad pixel-corrected and flat-fielded. Background sub- 
traction was made for each cube using the closest sky images. 
The images were rencentred using a Gaussian fit of the attenu- 
ated central star. The data were selected within each data cube 
using criteria based on the attenuated central star flux and the 
encircled energy between 0.4" and 0.55". The images were then 
binned every 6 s and derotated into a reference frame where 
the pupil was stabilized in order to simulate angular differential 
imaging (ADI, Marois et al. 2006). In this reference frame, the 
total field rotation provided by the manual offsets plus the natural 
pupi/field rotation is 17°. This data cube was then reduced using 
principal components analysis (PCA, Soummer et al. 2012), re- 
taining four components out of 105. 

The noise in the final reduced image was estimated using a 
sliding nine pixel-wide box to obtain a preliminary map of de- 
tection limits in magnitude. We corrected this map by computing 
the flux losses due to the PCA reduction. They were estimated by 
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Fig. C.l. Map of the detection limits in Jupiter masses set by the COND 
evolutionary models. The contours range from 60 to 150 M Jup with a 
step of 10 Mj up . 


injecting fake planets into the data cube at a 10-cr level and pro- 
cessing the data again. Last, these detection limits in magnitude 
were converted into detection limits in masses, using the COND 
(Baraffe et al. 2003) or BT-settl models (Allard et al. 2011), as- 
suming an age of 2 Gyr. The 2D-detection limits derived with 
the COND evolutionary models are presented in Fig. C. 1 . 
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